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Abstract. We show that a flame tracking/capturing scheme originally developed for deflagration fronts can be used to model 
thermonuclear detonations in multidimensional explosion simulations of type la supemovae. After testing the accuracy of the 
front model, we present a set of two-dimensional simulations of delayed detonations with a physically motivated off-center 
deflagration-detonation-transition point. Furthermore, we demonstrate the ability of the front model to reproduce the full range 
of possible interactions of the detonation with clumps of burned material. This feature is crucial for assessing the viability of 
the delayed detonation scenario. 
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1. Introduction 

The widely accepted standard model for type la supernovae 
(SNe la) is the thermonuclear explosion of a C-i-O white dwarf 
that has reached the Chandrasekhar m ass by means of mass ac¬ 
cretio n from a binary companion te.g. lHillebrandt & Niemev^ 
l200d and references therein). This scenario has recently re¬ 
ceived new support from the tent ative discovery of the compa n- 
ion star of Tycho’s supernova bv lRuiz-Lanuente Iil il2nn4 . 

While this idea has been around since the original work by 
iHovle & Fowled ( Il960l) and countless CPU cycles have been 
spent in search of a realistic model, some crucial aspects of 
the explosion physics are still not understood. Apart from the 
scientific challenge of reproducing the inner workings of these 
fascinating astron omical events of which we already know s o 
many details (e.g. iMatheson et al.lEf)04t iRenetti et al jl2004l> . 
their potential as cosmological distance indicators has pro¬ 
vided us with yet another strong motivation to construct predic¬ 
tive, self-consistent explosion models. Two subjects take center 
stage in most current debates: (i) the conditions under which 
the dynamical phase of t he explosion ignites a nd (ii) the ther¬ 
monuclear flame forms jWooslev et alJEoo3) . and the ques¬ 
tions surrounding the possibility of a deflagration-detonation- 
transition (DDT) jNiemever & Wooslevll997HKhokhlov et al.l 
Il997l) during the late p hase of the explosion, giving rise to 
a “delayed detonation” dWooslev & Weaved 1 1994 iKhokhlo^ 
Il99lt IVamaokaet alJll992ll . The modeling and consequences 
of the latter are the topics of this paper. 

Detection of intermediate mass elements in the early spec¬ 
tra of SNe la is usually interpreted in terms of an initial phase of 
subsonic thermonuclear burning (deflagration) that allows the 
star to expand to lower densities. The thin deflagration front 
is hydrodynamically unstable, most importantly as a result of 


the Raleigh-Taylor (RT) instability, so it quickly becomes fully 
turbulent. Dilferent models for the unresolved turbulent flame 
structure have been employed in multidimensional simulations 
of exploding white dwarfs E lemgTCr^ Hillehrandtlll99.4 
iKhokhlo^ 1 1 99-4 iReinecke et ^ 1999al) : however, the global 
properties of the explosion appear to be robu st with respect 
to the particular choice of subgrid-scale model (IReinecke et alJ 
l2002tlGamezoetal.l2003l) . 

On the other hand, interpretations of the simulations with 
regard to whether or not a DDT must occur differ consider¬ 
ably. Pure deflagration models with the currently achievable 
grid resolution produce too little ^®Ni when compared with 
one-dimensional models that successfully fit the observational 
data. Moreover, they generically leave behind unburnt material 
near the stellar core, whose amount is strongly constrained by 
the lack of carbon features in late SN la spectra. It has been 
claimed that this is a problem intrinsic to all turbulent deflagra¬ 
tion models which can only be solved by a delayed detonation 
that burn s up all of the C-hO le ft in clumps of varying sizes near 
the core dOamezo et alJ2003) . 

Owing to th e fine-tuning req uired to achieve a DDT in un- 
confin ed media (jN i^meverl l 9991 strengthened by the recent re¬ 
sults of iBell et alJ2004l) and to the promising trend in highly re¬ 
solved turbulent deflagratio n models for burning p rogressively 
more material in the core dNiemever et alJl2003li . we believe 
that the jury is still out on the question of delayed detonations. 
Nevertheless, we will demonstrate in this work that the level-set 
method that has been succuessfull y employed to model d efla- 
gration fronts in SN la simulations dReinecke et alJl999bh can 
be modified in a straightforward way to represent unresolyed 
detonations. Its main adyantage is its complete control oyer the 
front yelocity at each point, allowing us to mimic the interac- 
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tion of the detonation with unresolved and resolved clumps of 
ashes in a realistic manner. One of our main conclusions is the 
sensitivity of the results with respect to the way these interac¬ 
tions are implemented, supporting our claim that an accurate 
detonation model is necessary. 

In previous attempt s to simulate delayed detonations 
in mu ltiple dimensions (lArnett & Livnel[T994t lOamezo et all 
bOOdj) . the detonation front was unresolved and no model for 
the propagation speed was employed. While it is well-known 
that the velocity of planar detonations is robust with regard to 
numerical resolution, the same is not true in the presence of 
obstacles or clumps of ash (clos ely related are conse quences 
of the cellular front structure, see lTimmes et aiJ2nnnh . In fact, 
fully resolved calculations show that C-hO detonations can 
hardly pass through burned material at all, whereas this behav¬ 
ior is easily hidden by coarse resolution (lMaieill2005l) . Here, 
our approach is to compare the two limiting cases of (i) a deto¬ 
nation that passes through burnt material as a pure shock wave 
at the speed of sound and re-ignites on the other side and (ii) 
one that stalls immediately when it hits a clump of ash and 
needs to wrap around it to reach the other side. 

In all of our simulations, we let the code decide where 
and when the DDT should take place. The criterion for a 
DDT follows from the comparison of the thermal flame thick- 
ness with the Gibson leng t h of the turbulent flame b rush 
llNiemever & Wooslevlll99'^ iNiemever & Kersteinlll997h . As 

expected, we hnd that all DDTs occur near the outermost parts 
of the turbulent deflagration front, generally taking place ear¬ 
lier in off-center ignition models than in centrally ignited ones. 
It follows immediately that in order to burn the clumps of C-i-O 
near the core the detonation needs to propagate back into the 
star, having to pass through a foam of fuel and ashes along 
the way. This is the reason why an accurate treatment of the 
detonation-clump interactions is indispensable. 

In what follows we give a short description of the flame al¬ 
gorithm employed to model both the deflagration, as well as the 
detonation front. In Section|3 we discuss the results of a series 
of test calculations determining the suitability of our approach 
to represent supersonic burning fronts. In Section^ we present 
our DDT model and discuss the results and implications of our 
numerical study. Concluding remarks are made in Section im 

2. Level set deflagration and detonation model 

The numerical code applied in this work is base d on Reinecke’s 
version o f PROMETHEUS as described in iReinecke et al.l 
lll999albh . In partic ular, we empl oy the ley el set captur¬ 
ing/tracking scheme llSethianlll996l) used in these references 
for turbulent deflagration fronts to model an additional delayed 
detonation. In this approach, the thermonuclear burning front is 
associated with the zero leyel set of a linear distance function 
G(x, t). It thus represents an- 1-dimensional moying hypersur¬ 
face in an n-dimensional simulation. Eor a similar app r oach t o 
model deflagrations and detonations, see lEedkiw et alJ lll999l) . 

The front propagation is giyen by the temporal eyolution of 
G, 

— = (Vu-n + 5)|VG|, (1) 


depending on the fluid motion Vu of the unburnt state and 
the propagation speed s relatiye to it. The choice of s re¬ 
flects the kind of combustion front, whether deflagration (sub¬ 
sonic) or detonation (supersonic), and is proyided as an exter¬ 
nal parameter as described in the next paragraph. Details about 
the leyel set function a nd its implementation can be found in 
IReinecke et alJ lll999bl) . Eor the modeling of both deflagration 
and detonation fronts we used the cell ayeraged fluid yelocity v 
instead of Owing to this so-called pass iye implementation, 
as opposed to the full i mplern entation. cf. ISmilianowski et alJ 
lll997l) : lReinecke et alJlll999bl) . the transition between fuel and 
ashes is smeared out oyer about three grid cells. 

As in preyious simulations (IReinecke et alJfl999al l2002t) . 
the speed of the deflagration front is defined as the upper limit 
of the laminar and turbulent burning yelocities in the respectiye 
grid cell, idef = max(iiani, ■Stur)- Here, .iiam is a known function 
of co mposition, density, and temperature dTimmes & Woosle^ 
mu, while Stur is calculated by an additional subgrid-scale 
model taking fluid dynamics into a ccount on scales smaller 
than the computational grid size iNiemeyer & Hillebrandt 
|l22i- The detonation yelocities, sjet, are taken from ISharpe 
lll999t) as a function of the unburned fuel density. 

Behind the front all material is immediately transformed 
into the reaction products, whereas inside the mixed cells we 
determine the yolume fraction occupied by the burnt material 
and adapt the mass fraction of ashes correspondingly. Eor den¬ 
sities p > 5.25 X lO^g/cm^, the fuel is conyerted into ^®Ni, 
releasing a fusion energy of 7.86 x lO'^erg/g. At lower den¬ 
sities, the burning process only produces intermediate mass el¬ 
ements represented in our simulations by ^'^Mg and an energy 
release of 4.18 x 10'^ erg/g. Below p - 10’ g/cm^, no change 
of composition takes place. 


3. Detonation test calculations 

Since the front algorithm had so far only been applied for mod¬ 
eling deflagrations, its suitability for describing detonations 
had to be yerihed. Seyeral one- and two-dimensional test calcu¬ 
lations were carried out to inyestigate the influence of external 
factors like fuel density, grid resolution and time step on the 
front propagation. 

In all simulations, the unburnt material was characterized 
by a uniform composition with equal mass fractions of and 

and an initial temperature of T = 10^ K. Unless otherwise 
stated, all fronts were set to propagate in a positiye .r-direction 
(and a positiye {/-direction in the 2D case) on a Cartesian grid 
with a cell size of A = 10® cm. 


3.1. One-dimensional fronts 

We inyestigated planar fronts propagating on a 200 x 4 grid 
with reflecting boundaries at the left, top, and bottom edges 
of the computational domain and an outflow boundary to the 
right. Eigure^shows the prohles of temperature, density, fluid 
yelocity, and chemical composition of a detonation front 0.1s 
after its initialization. 

In order to inyestigate the ability of the front model to re¬ 
produce a giyen propagation yelocity idet, we compared the 
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Fig. 1. Temperature, density, fluid velocity, and chemical composition of a planar detonation front after 0.1 s. The spike in the 
density profile behind the shock is the result of a small numerical oscillation. 


time dependent x-positions of the detonation front with its the¬ 
oretically predicted behavior for three different initial densities 
p = 10^, 10*, and 10®g/cm*. To determine the position of the 
front, the values of the level set function in the direct neighbor¬ 
hood of the front were interpolated. 

Additionally, to verify the obtained x-values, we directly 
calculated the propagation velocity Sdet by averaging the mass 
fraction of ashes per cell over all 200 grid zones in the x- 
direction in every time step, multiplying with the physical grid 
length and dividing by time: 

Sdet = 200Xash(f)y. (2) 

Here, 2fash(0 denotes the averaged time dependent mass frac¬ 
tion of ashes. Both methods gave virtually identical results. For 
all three densities, the measured numerical velocity of the det¬ 
onation front idet exceeded the provided one, Sdet. by approxi¬ 
mately 10 %. As shown by further numerical tests, this behav¬ 
ior is affected neither by a change of the grid resolution nor by 
a considerable refinement of the time step, so that the error can 
be attributed to the inaccuracies in the passive implementation 
of the front model. 

Reconsidering the velocity profile in Fig. ^ the reason for 
the excessive propagation speed becomes clear. In our simula¬ 


tion, the fuel is at rest, while the ashes behind the front proceed 
with Vb ~ 2x 10* cm s“*. Averaging the two velocities results 
in a non-zero fluid velocity and therefore causes an error in Eq. 
njwhen replacing Vu by v. 

Since in these numerical tests the fuel had initially been 
at rest, the results shown above are only valid in the frame of 
reference of a non-moving gas. In order to test for a possible 
grid dependence of the front propagation scheme, the system 
of detonation front and gas was Galilei-transformed in a further 
simulation. Here, the C-i-O gas flowed uniformly in the negative 
x-direction with a velocity Ux - -1.16 x 10® cm s“*, and the 
reflecting boundary condition at x = 0 was changed to outflow. 
The front was ignited at the center of the grid. Since was 
chosen to be equal to Sdet at the density p - 10* g cm“*, the 
front should remain at rest in the grid frame on its right hand 
side while propagating with double velocity to the left. 

The profile of the front 0.018 s after its initialization can be 
seen in Fig.|2] The measured deviation (=s 7 %) is equal in both 
directions, showing that the front model is insensitive to Galilei 
transformations. This fact allows a simple correction of idet to 
account for the numerical error. 
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Fig. 2. Detonation profile 0.018 s after initialization in a C+O gas of velocity Vx = -1.16 x 10^ cm s 


3.2. Two-dimensional fronts 

For the investigation of fronts in two dimensions, the grid was 
expanded from 4 to 200 cells in the ^-direction. Calculations 
with circular detonations showed that the front kept its shape 
throughout the simulation time, differing insignificantly from 
the exact circle geometry (Fig.|3. As in the case of planar det¬ 
onations, the small deviations of the front velocity from its the¬ 
oretical value can be explained by the fact that the front is not 
advected with the speed of the unburnt matter in Eq. o but by 
the average speed in the burning cells. 

Com paring our findings to the results of iReinecke et al.l 
jl999bl) . one can see that the front algorithm causes errors like 
those for both kinds of burning fronts, with the difference that 
the speed is too low in the case of deflagrations while being 
too high in the case of detonations. The reason for this op¬ 
posite behavior is that in the earlier work, all material behind 
the front is at rest while the fuel expands at high velocities. 
Underestimation of the propagation speed is therefore a conse¬ 
quence of averaging the fluid velocities in the cells cut by the 
front. In contrast to our previous expectations, the higher den¬ 
sity jump across the front does not result in a larger inaccuracy 
for detonations. 


4. Two-dimensional simulations of delayed 
detonations 

The supernova simulations were carried out in cylindrical co¬ 
ordinates on a Cartesian grid in the r-z-plane. The center of the 
white dwarf was placed into the origin of a 256 x 256 mesh 
with equally spaced cells of a side width A = 10® cm. Thus, 
we modeled one quadrant of the star, assuming mirror symme¬ 
try along the r = 0 and z - 0 axes. To take the expansion of 
the star into account during the explosion and to avoid violat¬ 
ing of mass conservation due to an outflow of stellar material 
over the grid borders, the size of the outermost grid zones grew 
exponentially beginning with the 222nd zone. 

The white dwarf, constructed in hydrostatic equilibrium for 
a realistic equation of state, had a central density of 2 x 10® 
g cm“^, a radius of about 1.9 x 10* cm, and a mass of 2.8 x 
10** g. The initial mass fractions of C and O were chosen to 
be Xc,o = 0.5, and the total binding energy turned out to be 
-5.2 X 10*° erg. 

The deflagration phase of the models Z l, Z3, Bl, and 
B5 was initialized and evolved identically to IReinecke et al.l 
(Il999al) . For modeling the delayed detonation, a second level 
set function (LSet2) was implemented in addition to the first 
one (LSetl) describing the deflagration front. In order not to 






























I. Golombek and J.C. Niemeyer: Multidimensional Delayed Detonations 


5 



x[cm] 

Fig. 3. Front geometry of a centrally initialized circular front 
propagating outwards. The dashed line represents an exact cir¬ 
cle. 

disturb the flame propagation during the deflagration phase, the 
value of LSet2 was set to a negative value. 

A DDT was assumed to occur as soon as the transi¬ 
tion from flamelet to distributed burni ng regime took place 
for the first time during the expl osion llNiemever & WooslevI 
Il997t^emever & Kersteinll997h . In order to identify the cor¬ 
responding time and grid location, the thermal flame width 
d was compared to the Gibson length, defined by /cibs = 

A(4m/2^) , in each time step. The valu es of the flame width 

as a fu nction of density were taken from iTimmes & WooslevI 
lll992h . A detonation was initialized as a circular front with a 
radius of 10^ cm at the point where the condition /cibs = 5 was 
hrst satished. 

So far, the question whether a thermonuclear detonation 
front can propagate through burnt stellar material has not re¬ 
ceived much attention. It is conceivable that the pressure wave 
is strong enough to ignite the carbon again when it re-enters 
a region of fuel. However, fully resolved simulations suggest 
that this is not th e case even for relatively small clumps of 
ash (E aie ^ |2005|) . Depending on the behavior of the detona¬ 
tion front, isolated bubbles of fuel may be left inside the stellar 
core even after a delayed detonation has taken place, in dis¬ 
agreement with spectral observations in the nebular phase. 

In order to test the capability of our detonation model to 
cover the whole range of expected front behavior, we investi¬ 
gated the limiting cases of a detonation that first crosses the 
burnt material as a shock wave with sound velocity and then 
re-ignites on the other side (Case a), and a detonation that dies 
immediately upon running into ashes (Case b). The latter was 
modeled by setting = 0 in burnt regions. 

The time evolution of the front geometry and flow field in 
the non-central one-bubble model B1 is shown in Fig. 0]where 
we compare three snapshots of the explosion scenarios a and 
b at f = 0.7 s, 0.74 s and 0.84 s. As an additional measure for 


Table 1. Global characteristics of Case "a" models (detonation 
re-ignited after passing through clumps of ashes). 


model 

F 

^nuc 

[10^* erg] 

Ekin 

[10^^ erg] 

MMg 

[Mol 

Mni 

[Mo] 

Zla 

1.53 

0.90 

0.66 

0.61 

Z3a 

1.44 

0.87 

0.58 

0.65 

Bla 

1.75 

1.19 

0.87 

0.42 

B5a 

1.51 

0.95 

0.68 

0.53 


the expansion of the star, the central density is noted in every 
plot. By construction the deflagration phase is identical for both 
Cases a and b. 

After the DDT which takes place at f = 0.68 s at a dis¬ 
tance r ^ 1.5 X 10* cm, the differences between Models Bla 
and Bib become visible. The detonation in Model Bib cannot 
cross the bubble of burnt material until the latter breaks apart at 
its thinnest point at f ss 0.84 s (lower right plot of Fig. 0}. Notice 
that the breaking is only a numerical artifact since burnt regions 
cannot disconnect in the absence of flame quenching; more re¬ 
alistic, i.e. highly resolved and three-dimensional, simulations 
are needed to study the effects of ash clumps on predictions of 
the delayed detonation scenario. 

For completeness, the global properties of Models Zla to 
B5a are listed in Table [Q In all of these models, no isolated 
regions of fuel developed at late times, hence the Case b simu¬ 
lations were nearly identical to Case a in terms of global quan¬ 
tities. 

As can be seen from the energies and nickel amounts, all 
four models produce healthy explosions with high total ener¬ 
gies. The off-center ignition Scenarios B5 and especially B1 
are more powerful than the central explosion Models Z1 and 
Z3. One reason for this is that an off-center deflagration phase 
leaves large amounts of unburnt material at high densities near 
the center of the white dwarf, whereas a centrally ignited flame 
causes most material to expand significantly before the deto¬ 
nation sets in. In addition, the transition from deflagration to 
detonation took place earlier (0.6 s -0.7 s) in the off-center than 
in the central (~ 0.85 s) models. 

4.1. Conclusions 

In order to simulate delayed detonations in multidimensional 
explosion models of type la supernovae, we modified and 
tested a front algorithm previously developed to model sub¬ 
sonic deflagrations. The test calculations with detonation fronts 
presented in this paper show that, in spite of a signihcantly 
higher burning speed and density jump across the front, the 
inaccuracies in reproducing the correct front velocity are the 
same as those for deflagrations (< 10 %). In particular, the de¬ 
viations were shown to be independent of the flow velocity in 
the grid frame and therefore easily correctable by renormaliz¬ 
ing the prescribed front velocity. 

In a first exploratory set of two-dimensional explosion sim- 
ulations whose deflagr ation phases were identical to those in 
feeinecke et al.l ill999al) . we used two independent implemen¬ 
tations of the front model, one describing the turbulent flame 
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Fig. 4. Temporal evolution of the front geometry and velocity field for the Models Bla (top row) and Bib (bottom row). The 
snapshots were taken at equal times. The convoluted front corresponds to the initial deflagration, and the smooth curve that grows 
with time depicts the detonation front. 


front in the deflagration stage and the other accounting for 
the detonation front after an assumed deflagration-detonation- 
transition (DDT). The time and location of the DDT were de¬ 
termined by tracking the ratio of Gibson length and thermal 
flame thickness, and triggering the det onation where this pa- 
rametp dropped below unity (foll owing iNiemever & WooslevI 
ll997tlMemever & Kersteinll997ll . This procedure gave rise to 
DDTs at radii roughly 1.5 x 10® cm off-center and times of 
0.6 .. .0.8 s after ignition of the deflagration. Off-center defla¬ 
gration models generically resulted in higher energy release as 
a consequence of the earlier DDT and correspondingly higher 
core fuel density during the detonation phase. 

More realistic three-dimensional simulations are planned 
for future work. We expect a noticeable difference as a result 
of the changed flame front topology. However, the majority of 
scales of unburned clumps won’t be resolved in these simula¬ 
tions and will have to be modeled on subgrid scales. 

We emphasize the importance of using a tracking/capturing 
scheme for the detonation, as well as for the deflagration, as op¬ 
posed to the traditional approach to model detonations as un¬ 
resolved supersonic burning fronts. In addition to being able to 
use accurate tabulated detonation velocities as a function of the 
unburned fuel density, it allows full control over the interaction 
of the detonation with regions of burned material. The viabil¬ 
ity of the delayed detonation scenario hinges on the capability 


of the detonation to burn the remaining C-hO in the core, and 
hence on its ability to cross the foam of fuel and ashes in the 
Rayleigh-Taylor mixing region. These questions will be inves¬ 
tigated in a forthcoming publication. 
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